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Abstract 

In  this  paper  we  consider  models  for  noncausal  processes  consisting  of 
discrete-time  descriptor  dynamics  and  boundary  conditions  on  the  values  of  the 
process  at  the  two  ends  of  the  interval  on  which  the  process  is  defined.  We 
discuss  the  general  solution  and  well-posedness  of  systems  of  this  type  and 
then  apply  the  method  of  complementary  processes  to  obtain  a  specification  of 
the  optimal  smoother  in  terms  of  a  boundary-value  descriptor  Hamiltonian 
system.  We  then  study  the  implementation  of  the  optimal  smoother.  Motivated 
by  the  Hamiltonian  diagonal izat ion  results  for  non-descriptor  systems,  we  show 
how  the  descriptor  Hamiltonian  dyneimics  can  be  transformed  to  two  lower-order 
systems  by  the  use  of  transformation  matrices  involving  the  solution  of  two 
generalized  Riccati  equations.  We  present  several  examples  illustrating  our 
results  and  the  nature  of  the  smoothing  solution  and  also  present  equations 
for  covariance  analysis  of  boundary-value  descriptor  processes  including  the 
smoothing  error.  In  addition  we  discuss  several  open  problems  and  connections 
with  other  related  results. 
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I.  Introduction 

The  class  of  descriptor  systems  was  introduced  by  Luenberger  [1]  to 
describe  the  dynamics  of  certain  linear  systems  for  which  standard  state  space 
representations  are  not  particularly  natural  or  appropriate.  Since  their 
introduction  numerous  studies  have  been  performed  to  investigate  the 
properties  of  these  systems  and  the  solution  of  control  problems  for  them 
(see,  for  example,  [2]  -  [9],  [20],  [21]  and  the  references  cited  therein). 

The  fundamental  property  that  all  of  these  studies  have  had  to  deal  with,  in 
some  form  or  another,  is  the  fact  that  the  system  function  matrix  for  such  a 
system  is  not  proper,  leading  to  impulsive  behavior  in  continuous- time  and 
giving  rise  to  noncausal  responses  in  discrete-time.  The  noncausality  of 
these  models  makes  them  a  natural  choice  for  modeling  spatially  (rather  than 
temporally-) varying  phenomena,  and  in  this  context  it  is  natural  to  consider 
descriptor  models  with  general  boundary  conditions  rather  than  with  initial 
conditions  or  the  special  constrained  forms  for  boundary  conditions  found  in 
the  literature.  Indeed,  if  one  considers  generalizations  of  descriptor  models 
to  more  than  one  independent  variable,  one  finds  that  these  models,  together 
with  appropriate  boundary  conditions,  arise  in  many  contexts  such  as  in 
describing  random  fields,  electromagnetic  problems,  gravitational  anomalies, 
etc. 

The  investigation  of  standard  (i.e.  not  descriptor)  boundary-value  models 
in  one  independent  continuous  variable  was  initiated  by  Krener  [12]  -  [14]  who 
has  investigated  many  of  their  fundamental  properties.  Adams,  et  al.  [10] 
developed  a  general  approach  to  estimation  for  boundary-value  models  and 
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applied  it  in  [11]  to  develop  efficient  estimation  algorithms  for  processes 
described  by  the  model  introduced  by  Krener.  In  this  paper  we  extend  our 
estimation  methodology  to  two-point  boundary-value  descriptor  systems 
(TPBVDS’s),  i.e.  discrete-time  descriptor  models  in  one  independent  variable 
and  with  general  boundary  conditions.  To  our  knowledge  this  represents  the 
first  study  of  descriptor  models  devoted  to  estimation,  and  as  we  will  see, 
our  analysis  uncovers  both  some  important  similarities  and  differences  with 
estimation  problems  for  standard  state  space  models  and  several  important 
problems  whose  solutions  remain  for  the  future.  These  questions  have  in  fact 
inspired  the  development  of  a  system-  theory  for  TPBVDS's  [25],  several 
elements  of  which  will  be  used  in  the  present  development.  Furthermore,  in 
another  paper  [15]  we  use  the  results  developed  here  in  our  investigation  of 
efficient  estimation  algorithms  for  random  fields  describable  in  terms  of  a 
particular  class  of  boundary-value  descriptor  systems  in  two- independent 
variables. 

In  the  next  section  we  introduce  the  class  of  TPBVDS’s  and  perform  some 
preliminary  analysis.  In  particular,  we  discuss  the  wel 1-posedness  of  such  a 
system  and  a  general  method  of  solution  for  TPBVDS’s.  In  Section  III  we  apply 
the  results  of  [10,  26]  to  the  fixed-interval  smoothing  problem  for  an 
nth-order  TPBVDS.  As  we  show,  aside  from  a  boxmdary  effect  which  can  be  dealt 
with  separately,  the  resulting  smoother  is  itself  naturally  described  as 
TPBVDS,  in  this  case  of  dimension  2n.  In  Section  IV  we  address  the  question 
of  implementation  of  the  smoother.  Motivated  by  the  "Hamiltonian 
diagonal ization"  results  in  [11,  22]  for  non-descriptor  systems,  we 
investigate  two  procedures  for  forward-backward  diagonalization  of  the 
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smoother  equations.  These  procedures,  which  are  illustrated  in  Section  V, 
point  out  connections  with  other  work  on  descriptor  systems  and  also  lead  to 
several  solved  and  open  problems  related  to  generalizations  of  causal 
system- theoretic  concepts  to  TPBYES’s.  These  are  presented  and  discussed  in 
Section  VII  following  our  analysis  of  the  smoothing  error  in  Section  VI. 
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II.  Two-Point  Boxindary-Value  Descriptor  Systems 

The  TPBVDS  considered  in  this  paper  satisfies  the  difference  equation 
Ex(k+1)  =  Ax(k)  +  Bu(k)  (2.1) 

with  the  two-point  boundary  condition 

VqX(O)  +  Vj^(K)  =  V  (2.2) 

Here  u(k)  is  eui  mxl  input  sequence  defined  on  the  discrete-time  interval 
[0,  k-1],  x(k)  is  the  n-dimensional  boundary  value  process,  v  is  the  n-vector 
of  boundary  values,  and  E,  A,  B,  V^,  and  are  matrices  of  appropriate 
dimensions.  Furthermore  we  assume  that  {E,  A)  form  a  regular  pencil  (i.e. 
|zE-Al  XO). 

As  in  [2],  we  can  rewrite  (2.1),  (2.2)  as  a  single  set  of  equations 
yk  =  26u  (2.3) 

where 

X-  =  (x*(0),...x'(K))  (2.4a) 

u'  =  (u'(0) . u’(K-l),v‘)  (2.4b) 

E  0  . 

A.  E  0  . 

•  • 

•  • 

. o'  -K 

0 . 0 


0 

0 


'K 


(2.5a) 


S  =  diag  (B . B,  I)  (2.5b) 

We  see  from  this  immediately  that  the  well-posedness  of  (2.1),  (2.2)  is 
equivalent  to  the  invertibi lity  of  iP.  Much  more  can  be  said  about 
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well-posedness  and  the  solution  of  (2.1).  (2.2),  and  we  refer  the  reader  to 
[25]  for  details.  We  limit  ourselves  here  to  describing  one  method  for 
solving  (2.1),  (2.2)  that  provides  us  with  an  alternate  well-posdness 
condition  and  with  a  method  for  the  implementation  of  the  smoother  developed 
in  Sections  III  and  IV. 

To  begin,  from  Kronecker’s  canonical  form  for  a  regular  pencil  [17]  we 
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can  find  nonsingular  matrices  T  and  F  so  that 


I  0 

°  V 


(2,6) 


fat"^  = 


0 

0  I 


(2.7) 


and  so  that  all  of  the  eigenvalues  of  A^  and  A^  have  magnitudes  no  larger  than 
1.  Furthermore  if  1zE-A|  has  no  zeros  on  the  unit  circle,  then  all  of  the 
eigenvalues  of  A^  and  A^  are  strictly  inside  the  unit  circle.  In  this  case  we 
will  say  that  {E,A}  is  f oryfard-backward  stable. 


uie  decomposition  in  [17]  splits  the  pencil  zE-A  into  forward  dynamics 
corresponding  to  a  pencil  of  the  form  zI-A^  and  backward  dynamics 

-I  ~  ~ 

corresponding  to  z  -A^^  where  Aj^  is  nilpotent.  The  only  difference  in  (2.6), 

(2.7)  is  that  the  unstable  forward  modes  of  A^  have  been  shifted  into  the 
backward  dynamics  A^. 


=  Tx(k) 


(2.8) 


x^(k) 


Then,  we  obtain 

Xj.(k+1)  =  A^x^(k)  +  B^u(k) 
=  -^j^Ck+l)  -  B^u(k) 

where 


=  FB 


(2.9a) 

(2.9b) 


(2.10) 


and  (2.9a),  (2.9b)  are  asymptotically  stable  recursions  if  (E,  A)  is 
forward-backward  stable.  Finally,  given  the  transformation  (2.8),  the 
boundary  condition  (2.2)  takes  the  form 


(^f,o:\,o]  " 


Xf(0) 

+  rv„ , 'v, 

Xf(K) 

>S(0)_ 

L  f ,K.  b,Kj 

_Xb(K) 

=  V 


(2.11) 


(2.12) 


Employing  the  forward/backward  representation  (2.9)  of  the  dynamics,  a 
general  solution  to  (2.1),  (2.2)  is  derived  as  follows.  Let  x^^(k)  denote  the 
solution  to  (2.9a)  with  zero  initial  condition,  and  let  x^*^(k)  denote  the 
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solution  of  (2.9b)  with  zero  final  condition.  Then 


x^(k)  =  A^^x^CO)  +  x^°(k) 


(2.13a) 


xj^(k)  =  \(K)  +  Xj^°(k) 


(2.13b) 


Substituting  (2.13)  into  (2.11)  and  solving  for  x^(0)  and  Xjj(K)  yields 


x^(0) 

x^(K) 


=  H-'{v  -  Vj  ,jKj°(K)  -  V^,oX^°(0)) 


(2.14) 


where 


(2.15) 


Finally,  substituting  (2.14)  into  (2.13)  we  obtain 


x^(k) 

A^^  0 

.-1 ,  0,  .  0,  . . 

On 

x^  (k) 

«  -  "^f.K^f  ^  -  \,0^  ^ 

On  ^ 

(2.16) 


The  solution  in  the  original  basis  can  then  be  obtained  by  inverting  (2.8). 

Asstuning  that  {E,  A)  is  forward-backward  stable,  the  solution  procedure 
is  just  described  consists  of  stable,  forward/backward  recursive  computations 
for  Xp^,  followed  by  the  correction  for  the  actual  boundary  conditions 
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given  by  the  first  term  on  the  right-hand  side  of  (2.16).  Note  also  that  this 
procedure  also  provides  us  with  another  necessary  and  sufficient  condition  for 
the  well-posedness  of  (2.1),  (2.2),  namely  the  invertibility  of  H  in  (2.15). 
This  condition  is  the  analog  of  that  described  by  Krener  [12]  -  [14]  for 
standard  boundary-value  problems.  Note  that,  as  one  would  expect,  not  all 
choices  of  boundary  conditions  lead  to  well-posed  problems,  and  the  conditions 
that  Vq  and  must  satisfy  depend  heavily  on  the  structure  of  E  and  A.  For 

example,  as  is  well  known,  the  initial  value  problem  (Vq  =  I,  Vj^  =  0)  is  not 
well-posed  if  E  is  singular.  This  can  easily  be  seen  from  (2.15)  or  from 
(2.5a),  since  the  last  block  of  columns  then  is  not  of  full  rank. 
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III.  The  Optimal  Smoother 

Consider  now  a  stochastic  process  x{k)  satisfying  (2.1),  (2.2)  (which  we 
assume  is  well  posed)  where  u(k)  and  v  are  independent,  zero  mean  and 
Gaussian,  v  has  covariance  H^,  and  u(k)  is  a  white  sequence  with  covariance  Q. 
In  this  section  we  examine  the  estimation  of  x(k)  given  the  interior 
observations 

y(k)  =  Cx(k)  +  r(k),  k  e  [1,  K-1]  (3.1) 

and  the  boundary  measurements 

y^  =  WqX(O)  +  Wj^(K)  +  r^  (3.2) 

Here  r(k),  r^,  u(5),  and  v  are  mutually  independent,  r^  is  zero  mean  Gaussian 
with  covariance  IT^,  and  r(k)  is  zero  mean,  Gaussian,  and  white  with  covariance 
R. 


In  order  to  derive  the  optimal  smoother,  we  introduce  notation  analogous 
to  (2.4),  (2.5) 

y  =  ‘foe  +  r  (3.3) 

where 


y'  =  [y'(i).  y'(2) . y’(K-i),  y^-] 

r’  =  [r'(l),  r'(2),...,r'(K-l),  r^*] 


0 

0 

0 


(3.4a) 

(3.4b) 


(3.5) 
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Also,  the  covariances  of  u  in  (2.4b)  and  r  in  (3.4b)  are  given  by 

Q  =  diag(Q . Q.  U^)  (3.6a) 

^  =  diag  (R _ .R.  H^)  (3.6b) 

Our  problem,  then  is  to  estimate  x  given  y,  and  the  approach  we  adopt  is 
the  method  of  complementary  processes  Introduced  in  [26]  and  elaborated  upon 
in  [10,  11].  Specifically,  suppose  that  we  can  construct  a  random  vector  z 
that  is  complementary  to  y  in  the  sense  that  (i)  it  is  independent  of  y  and 
(ii)  the  transformation  from  (u,r)  to  (y.z)  is  linear  and  invertible.  Then  we 
can  write  x  explicitly  as  a  linear  function  of  y  and  z,  and,  thanks  to  (i)  can 
obtain  x  simply  by  setting  z  to  zero.  In  the  present  context,  since  x  is 
specified'  implicitly  by  (2.3),  we  also  obtain  and  implicit  representation  for 
z.  Specifically,  as  we  verify  below,  z  is  given  by  the  following 

(3.7) 

z  =  -S'X  +  C  ■"u  (3.8) 

where 

X’  =  [\'(1).....X’(K).X'(0)]  (3.8) 

(the  reason  for  our  particular  choice  of  labeling  of  components  in  (3.8)  will 
be  made  clear  shortly).  Note  that  (3.7)  also  has  an  interpretation  as  a 
TPBVDS,  but  we  defer  discussion  of  this  until  our  related  discussion  of  the 
smoother  itself. 

As  a  first  step  in  verifying  (3.7),  (3.8)  note  that  (3.7)  is  well-posed 
since  is  invertible.  Next  note  that  the  independence  of  y  and  z  can  be 
obtained  by  direct  computation: 

E{yz’}  =  E{['ef“^ah  +  r][-£§'(if)“^*2S~^r  +  C~^u]'} 

=  -  ‘ey’”^26  =  o 


y”X  =  <€'2^  ^r 


.-1 


(3.9) 
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We  next  show  that  we  can  compute  x  and  X  from  y  and  z.  Specifically,  using 
(3.8)  to  eliminate  u  and  (3.3)  to  eliminate  r,  we  find  that  (2.1),  (3.7)  are 
equivalent  to 


•  • 

X 

'mz 

X. 

(3.10) 


The  matrix  on  the  left-hand  side  of  (3.10)  can  be  shown  to  be  invertible  as 
follows.  Since  if  is  invertible,  we  need  only  show  that  the  Schur  complement 

D  =  y  +  (3.11) 

is  invertible.  Note  that 

D(y”)"^  =  I  +  MK  (3.12) 

where  M  =  ‘6’^  >  0  and  K  =  y  ^^S§’(y”)  0.  The  invertibility  of  D  then 

5 

follows  from  the  fact  that  MK  cannot  have  negative  eigenvalues.  Filially,  once 
we  have  recovered  x  and  X  from  y  and  z,  u  and  r  can  be  obtained  from  (3.8)  and 
(3.3),  respectively. 

Next,  by  setting  z  to  zero  in  (3.10)  we  obtain  the  implicit  equations 
defining  the  optimal  smoothed  estimate  x: 


if 


'e  if' 


X 

0 

X, 

(3.13) 
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Suppose  MKv  =  Xv. 


Then  v'K'MKv  =  Xv'K’v,  so  that  X  =  (v'K’MKv*/(v'K'v)  >  0. 


13. 


This  again  defines  a  well-posed  TPBVDS,  but  to  obtain  the  most  illuminating 
form'  of  this  system  requires  a  permutation  of  the  equations  and  variables  in 
(3.13).  Specifically,  it  is  straightforward  to  verify  that  (3.13)  is 
equivalent  to 


Sff  =  T7 


where 

f  =C(x'(0).  V(0)).  (x'(l).  X'(l)).....(x'(N).X'(N))] 


(3.14) 


(3.15a) 


0 


\c‘R  V(N-i; 


(3.15b) 


0  -o . 0  V 

S  0 . .0  0 

-d  5  . 0  0 


0  0  0  0  -d  8 

■^21  °  °  . ^'^22 


(3.16) 


with 


E 

-BOB'" 

A 

0  ■ 

s  = 

.  d  = 

—  1 

.0 

-A' 

.-C’R  C 

-E'. 

< 

1 

> 

o 

_ ! 

0 

0 

-1 

t  = 

’  12 

-1 

%  \  Vo'J 

ro^b  \ 

0 

'^0 

-n  1 

V 

0  ■ 

V  ' 

K  J 

1f  = 

’  22 

E' 

(3.17) 


(3.18a) 


(3.18b) 
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Comparing  the  form  of  ?  in  (3.16)  to  that  of  if  in  (2.3).  we  see  that 

(3.14)  is  almost  a  standard  TPBVDS  except  for  the  top  row  of  equations  -  i.e. 

the  fact  that  11^^  in  (3.16)  appears  rather  than  -d  and  that  '^■^2  present  at 

all.  This  is  a  consequence  of  the  discrete  nature  of  the  time  index  and  the 

0 

intrinsic  asymmetry  of  the  model  (2.1),  (2.2).  We  can,  however,  reduce  these 
equations  to  a  standard  TPBVDS  by  means  of  a  basic  technique  in  the  analysis 
of  boundary-value  systems  [14,25].  Specifically,  we  can  think  of  (3.14)  as  a 
TPBVDS  with  boundary  values  consisting  of  (x'(0).  X'(0))'  and  (x'(N),  X'(N))‘. 
Because  of  the  well-posedness  of  (3.14)  it  is  possible  to  eliminate  some  of 
the  variables  from  (3.14)  by  solving  for  them  in  terms  of  the  remaining 
variables.  More  specifically,  it  is  possible  to  move  the  boundary  values 
inward  by  eliminating  boundary  values  at  one  end  of  the  interval,  the  other, 
or  both.  One  can  iterate  this  process,  and  in  fact  this  type  of  recursion 
forms  the  basis  for  a  notion  of  state  for  boundary  value  systems  [14,25].  For 
our  purposes  here,  however,  we  need  only  consider  a  single  step  of  this  type. 
Specifically,  the  invertibility  of  2?  implies  that 


t 


11 


t, 


21 


0 

Note  that  u(k)  is  defined  on  [0,  K-1],  while  x(k)  is  defined  on  [0,K]. 
Referring  to  [10],  it  is  not  possible  in  the  discrete  index  case  to  define  the 
domain  on  which  x  and  u  are  defined  and  the  boundary  of  that  domain  so  that 
either  the  boundary  is  contained  in  or  disjoint  from  the  domain. 
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has  full  column  rank  and  thus  that  we  can  eliminate  (x’ (0) ,X(0) ) ’  as  follows. 
We  construct  matrices  and  such  that  [Mj,M23  has  full  row  rank  and 


[Ml-  M2] 


t 


11 


t, 


21 


=  0 


(3.19) 


If  we  then  premultiply  (3.14)  by  the  following  full-rank  matrix 


0  10. 

0  0  I  . 

0  0  0. 

Mj  0  0  . 


0  0 
0  0 

I  0 
0  M2 


we  obtain  a  TPBVDS  of  a  form  exactly  as  in  (2.1),  (2.2).  Specifically,  this 
computation  yields 


x(k+l) 

=  d 

x(k) 

+ 

X(k+1) 

•  . 

X(k) 

0 

<  •  "o”!  . 


k=l . K-1 


(3.20) 


with  bovindary  conditions 


Mjg 


x(l)' 

.Ml). 


t“l''l2 


x(N)- 

.MN). 


0 

“  1 

=  Mi 

•  ■ 

"“2 

* - 

(3.21) 
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By  construction  we  know  that  this  system  is  well-posed.  Also,  once  we  have 

computed  x(k) ,  X(k) .  k=l . K,  we  can  determine  the  previously  eliminated 

boundary  values  x(0) ,  X(0) : 


x(0)- 

i 

0 

0 

.MO). 

■"o'V'j’b- 

+  1! 

21 

-If  '  S 
^11  ^ 


■x(l) 

LMi) 


~  1*^12  ^ 


21  22- 


x(K)l  ) 

.MK)J  j 


(3.22a) 


where 


D  =  +  ”<^21 ’*^21^’^  (3.22b) 

As  a  final  comment,  we  note  that  on  examination  of  (3.20),  (3.21)  and  the 
form  of  S  and  d  in  (3.17),  we  see  that  what  we  have  derived  is  a 
generalization  of  the  Hamiltonian  form  of  the  optimal  smoother  for  causal 
systems  (see,  e.g.  [11,22]).  This  immediately  suggests  the  possibility  of 
generalizing  methods  for  solving  smoothing  equations  such  as  diagonal ization 
of  the  Hamiltonian  dynamics  [11,22]  to  produce  forward  and  backward 
recursions.  Such  an  approach  is  described  in  the  next  section. 
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IV.  Implementation  of  the  Smoother 

In  this  section  we  discuss  several  appraoches  to  solving  the  smoothing 
TPBVDS  (3.20),  (3.21).  One  obvious  method  of  solution  is  the  direct 
application  of  the  method  described  in  Section  2  for  solving  general  TPBVDS’ s. 
The  question  that  then  arises  is  the  construction  of  the  similarity- 
transformations  that  block-diagonal ize  i  and  d  as  in  (2.6),  (2.7).  One 
obvious  answer  to  this  is  to  use  the  general  procedure  in  [17]  for  the 
computation  of  the  Kronecher  form  of  [S,d).  A  second  is  to  consider 
generalizations  of  Hamiltonian  diagonal ization  procedures,  which  are  developed 
in  the  following  two  subsections.  In  the  first  of  these  we  closely  parallel 
the  approach  used  for  non-descriptor  systems  and  are  led  to  descriptor  Riccati 
equations  and  decoupled  descriptor  dynamics.  As  we  will  see,  this  approach 
does  not  always  work,  and  this  leads  us  to  a  slightly  different  approach  in 
Section  4.2  involving  a  different  type  of  generalized  Riccati  equation  and 
producing  decoupled  non-descriptor  dynamics.  Open  questions  remain  concerning 
existence  of  solutions  to  these  equations,  but  as  we  discuss  in  this  and  in 
subsequent  sections,  this  approach  has  much  promise  emd  also  appears  to  point 
the  way  to  developing  the  relationship  between  system  -  theoretic  concepts 
such  as  reachability  and  observability  and  properties  and  eigenstructure  of 
the  smoother. 

4.1  Hamiltonian  Diagonal ization:  Method  1 

The  general  concept  of  Hamiltonian  diagonal ization  is  as  follows.  We 
seek  two  sequences  of  matrices,  M(k)  and  N(k)  so  that 

Ef(k)  0 
0  A^(k) 


M(k)£N~^(k+l)  = 


(4.1) 
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and 


M(k)gN~^(k) 


\{k)  0 

0  Eb(k) 


(4.2) 


In  this  case  the  2n-dimensional  descriptor  dynamics  of  (3.20)  can  be  decoupled 
into  two  n-dimensional  descriptor  systems  (coupled,  of  course,  through  the 
boundary  cond i t i ons ) . 

The  choice  of  the  sequences  M(k)  and  N(k)  is  far  from  unique,  and  the 
general  algebraic  equations  that  the  nxn  blocks  of  M(k)  and  N(k)  must  satisfy 
are  presented  in  [18]  and  [19].  In  this  subsection  we  present  one  choice  that 
is  the  direct  counterpart  of  the  method  used  in  [11]  for  non-descriptor 
continuous-time  boundary  value  processes  and  that  involves  descriptor  Riccati 
eqtiations  that  have  appeared  elsewhere  in  the  literature.  Specifically, 
suppose  that  P(k)  and  9(k)  are  invertible  matrix  sequences  satisfying, 
respectively,  the  following  forward  and  backward  descriptor  Riccati 
recursions: 

EP(k+l)E'  =  A[P"^(k)  +  C‘R"^C]"^A’  +  BOB’  (4.3) 

E'0(k)E  =  A'[0“^(k+1)  +  BQB’]"^A  +  C’R’^C  (4.4) 

In  the  case  of  causal  systems  (with  E  =  I),  (4.3)  is  the  recursion  satisfied 
by  the  one-step  forward  prediction  error  variance,  while  (4.4)  is  the 
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recursion  satisfied  by  the  inverse  of  the  backward  filtered  error  variance. 
Also,  define 

Z{k)  =  E'0(k)E  +  P'^(k)  (4.5) 

In  the  causal  case  and  with  appropriate  choices  of  initial  condition  for  P(k) 
and  final  condition  for  0(k) .  Z(k)  is  the  inverse  of  the  smoothing  error 
variance. 

Def ine 


M{k) 


Z  ^(k)A‘[  0  ^(k+l)+BQB']  ^ 


A[P  ^{k)+C’R  ^C]"^ 


-  Z"^(k) 


(4.6) 


N"^(k) 


Z~^(k)  P(k)E‘ 

-  6(k)EZ"^(k)  •  I 


(4.7) 


Some  algebraic  manipulations  verify  that  M(k)  and  N(k)  are  invertible  if  P(k) , 
0(k) ,  and  Z(k)  are,  and  if  we  perform  the  computations  involved  in  (4.1), 

(4.2)  and  define 


f(k)1  A 

=  N(k) 

Jj(k). 


'x(k)' 

.Mk). 


(4.8) 


7  -1 

The  acttial  quantities  P(k)  and  0  (k)  have  these  interpretations  only  if  the 

initial  and  final  conditions  P(0)  and  0(K)  are  appropriately  chosen.  In  this 

case  [P  ^(k)  +  C'R  ^C]  ^  is  the  forward  filtered  error  covariance,  while 

0  ^(K)  +BQB'  is  the  one-step  backward  prediction  error  covariance. 
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the  smoother  dynamics  (3.20)  decouple  into 

EP(k+l)f(k+l)  =  ACP‘^{k)+C’R~^C]"^  [f(k)+C‘R“V(k)]  (4.9) 

P(k)E'Tj(k)  =  [P"^(k)+C'R"^C]“^A'T7(k+l)+Z"^(k)C-R"V(k)  (4.10) 

The  boundary  conditions  in  the  transformed  coordinates  can  be  determined  from 
(3.21),  (3.10),  and  (3.12). 


As  a  simpler  alternative  one  might  consider  constant  transformations  M 


and  N  as  in  (4.6)  and  (4.7),  but  using  solutions  to  the  steady-state 
descriptor  Riccati  equations 


EPE'  =  A[P"^+C‘R  +  BOB' 

E*  0  E  =  A’[  0  +  BQB']“^A  +  C*r“^C 


(4.11) 

(4.12) 


Note  that  in  this  case  the  transformed  smoother  dynamics 

EPf(k+l)  =  A[P"^+C’R~^C]"^  Cf(k)+C'R"V(k)]  (4.13a) 

PE’77(k)  =  [P“^+C’R"^C]"^A’T7(k+l)+Z“^C'R~V(k)  (4.13b) 

involves  two  pencils  {Ej^,Aj}  =  (EP,  A[P  +C'R  C]  }  and  {E2,A2}  = 

{PE',  [P  ^+C'R  ^C]  ^A’}  that  are  transposes  of  one  another.  In  this  case  if 
we  follow  the  solution  procedure  outlined  in  Section  2.1,  if  the  matrices 
and  transform  {Ej^.A^}  into  the  form  shown  in  (2.6),  (2.7),  then  ^2  “  ^1 ' 
and  ^2  =  '  do  the  same  for  {E2,A2}. 


The  descriptor  Riccati  equations  we  have  introduced  have  appeared  in  the 
literature.  In  the  case  in  which  E  is  nonsingular,  which  was  studied  by  Laub 
in  [24],  it  is  clear  that  these  are  no  difficulties  in  solving  (4-3),  (4.4)  or 
equivalent  versions  not  involving  inversions  of  P  and  0)  nor  in  obtaining 
controllability  and  observability  conditions  under  which  (4.11),  (4.12)  have 
unique  positive  definite  solutions.  Furthermore  in  this  case  it  is  also 
possible  to  parallel  the  approach  in  [11]  (for  the  non-descriptor  case)  in 
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choosing  boundary  conditions  P(0)  and  0(K)  for  (4.3),  (4.4)  so  that  the 

boundary  conditions  associated  with  (4.10)  are  minimally  coupled.  Similarly, 

in  the  case  in  which  A  is  invertible,  we  can  do  something  analogous,  leading 

to  a  pair  of  dual  Riccati  equations,  essentially  by  reversing  time  (k  K-k) 

thereby  interchanging  the  roles  of  A  and  E.  While  the  approach  outlined  in 

8 

this  section  (or  its  dual)  works  when  either  A  or  E  is  invertible,  the 
difficulty  arises  when  both  E  and  A  are  singular.  As  pointed  out  by  Bender, 
singularity  can  cause  equations  such  as  (4.3)  to  fail  to  have  solutions  for 
particular  initial  conditions.  Also,  as  we  illustrate  through  an  example  in 
the  next  section,  when  E  and  A  are  both  singular  (4.11),  (4.12)  have  solutions 
only  in  an  uninteresting  case.  What  is  therefore  required  is  a  different 
approach-  Previous  studies  of  control  problems  for  continuous  or  discrete 
descriptor  systems  [9],  [21],  [23]  have  circumvented  this  difficulty  by 
deriving  and  dealing  with  lower-order  standard  Riccati  equations  (of  dimension 
equal  to  the  rank  of  E) .  In  our  case,  however,  we  are  interested  in 
diagonalizing  the  Hamiltonian  dynamics.  As  we  develop  in  the  next  section, 
this  is  possible  if  we  introduce  equations  that  are  not  quite  standard  Riccati 
equations  but  are  far  closer  to  them  than  (4.11),  (4.12). 

4.2  Hamiltonian  Diagonal ization:  Method  2 

In  this  subsection  we  focus  completely  on  time-invarint  versions  of  the 
transformations  (4.1),  (4.2).  The  key  to  the  transformations  are  the 

g 

Note  that  5  and  d  are  both  singular  if  either  E  or  A  is,  so  that  the 
procedure  in  this  section  does  work  on  a  class  of  nontrivial  Hamiltonian 
descriptor  dynamics.  See  [27],  [28]  for  investigations  of  discrete-time 
algebraic  Riccati  equations  by  examination  of  the  pencil  defined  by  S  and  d 
when  E  =  I. 
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generalized  Riccati  equations 

0  =  A’(E0"^E*  +  BQB')“^A  +  C'R"^C  (4.14) 

^  =  A(E'^~h  +  C'R"^C)~^A‘  +  BOB'  (4.15) 

Note  that  these  equations  are  "almost"  standard  Riccati  equations,  except  for 
the  presence  of  E  and  E'  multiplying  0  ^  eind  'I'  ^  in  the  terms  in  parentheses. 
While  there  appears  to  be  some  asymmetry  in  the  roles  played  by  E  and  A,  this 
is  an  illusion,  as  can  be  seen  by  introducing  an  additional  pair  of  matrices, 
specifically,  if  we  define 


S  =  E0  ^E’  +  BOB' 

(4.16a) 

we  see  that 

0  =  A’S"^A  +  C‘R"^C 

(4.16b) 

Similarly,  by  introducing 

T  =  E’^^'^E  +  C‘R"^C 

{4.17a) 

we  obtain 

=  AT"^A’  +  BOB' 

(4.17b) 

Consequently,  we  can  view  (4.16)  and  (4.17)  individually  as  pairs  of  equations 
to  be  solving  for  (S,0)  and  (T.'l'),  respectively.  We  assume  throughout  this 
section  that  positive  definite  solutions  for  these  four  quantities  exist.  As 
in  the  previous  section,  if  either  E  or  A  is  invertible,  we  can  reduce  these 
equations  to  standard  Riccati  equations  and  therefore  cam  obtain  the  usual 
type  of  reachability  and  observability  conditions  for  existence  of  such 
solutions.  Also,  as  we  illustrate  in  the  next  section  these  equations  admit 
positive  definite  solutions  even  in  cases  in  which  both  E  and  A  are  singular. 
General  conditions  for  existence  and  uniqueness  of  positive  definite  solutions 
remain  open,  and  in  Section  VII  we  briefly  discuss  this  and  several  related 


questions. 
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Consider  next  the  matrices 


I  AT 

A'S“^  -I 


N  = 


E  -'P' 

9  E. 


(4.18a) 


(4.18b) 


The  invertibility  of  N  is  immediate  from  the  invertibility  of  0  and  the 
invert ibility  of  the  Schur  complement 
-'P  -  E*e"^E 

Similarly  the  invertibility  of  M  follows  from  the  invertibility  of  -I  and  of 
the  Schur  complement 
I  +  A’S'^AT 

(which  is  invertible  since  T  >  0.  A'S  ^A  >  0  so  that  the  eigenvalues  of 
A’S  ^AT  are  nonnegative). 

It  is  a  straightforward  exercise,  using  (4.16),  (4.17)  to  show  that 


MSN  ^  = 


0  A'S  ^E0  ^ 


Ms^N 


-1 


''AT  o' 


0 


(4. 19a) 


(4. 19b) 


Therefore,  if  we  premultiply  (3.20)  by  M  and  make  the  change  of  coordinates 


6(k)- 

x(k)' 

=  N 

A 

b(k)J 

LMk)J 

(4.20) 
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the  smoother  dynamics  are  transformed  into  s tandard  non-descriptor  recursions  •' 
6(k+l)  =  +  AT“^C'R"V(k)  (4.21a) 

^(k)  =  A’S"^E0"^^(k+l)  +  C*R"V(k)  (4.21b) 

with  boundary  conditions 


■5(1)- 

—  1 

■a(K)- 

.^(1). 

0 

+ 

■  • 

0 

2 

(4.22) 


Note  that  (4.21)  consists  of  a  forward  recursion  (a)  and  a  reverse  recursion 
(b) .  with  coupled  boundary  conditions  (4.22).  The  approach  outlined  in 
Section  II  (see  (2.13)  -  (2.16))  can  then  be  used  directly  to  obtain  the 

solution.  Once  this  is  accomplished,  we  cam  recover  x(k)  aind  X(k) ,  k=l . N 

by  inverting  (4.20),  i.e.  from  the  relationship 

x(k)  =  [9  +  E'\?''^E]"^[^(k)  +  E‘'J^"^6(k)]  (4.23a) 

X(k)  =  'j/"^E[0  +  E’'l>"^E]"^^(k)  -  +  E0~^E']"^5(k)  (4.23b) 

/s  yv 

and  then  can  recover  x(0) ,  X(0)  from  (3.22).  Note  that  since  one  is  generally 
interested  only  in  x.  it  is  only  necessary  to  solve  for  X(l)  and  X(N)  in 

A 

(4.23b)  in  order  to  be  able  to  determine  x(0)  from  (3.22). 
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V.  Examples 

In  this  section  we  first  present  an  example  illustrating  our  smoothing 
results  for  TPBVDS’s  and  then  introduce  the  class  of  cyclic  systems  in  a 
second  example. 

Example  5.1*  As  we  indicated  in  the  previous  section,  the  case  in  which 
either  E  or  A  is  invertible  can  be  thought  of  as  a  slight  generalization  of 
the  causal  case  (perhaps  with  time  reversal),  and  consequently  both  of  the 
Riccati-like  methods  of  the  previous  section  (or  the  dual  of  the  method  of 
Section  4.1)  work  without  difficulty.  In  this  example,  we  look  at  a  system 
for  which  both  E  and  A  are  singular  and  first  illustrate  the  problems  with  the 
method  of  Section  4.1  and  the  apparent  superiority  of  the  approach  in  Section 
4.2. 

Consider  the  descriptor  system  with 


'i  o' 

0  o' 

E  = 

A  = 

o 

.0  1. 

In  this  case  it  is  not  difficult  to  check  that  difficulties  arise  in  solving 
the  time-varying  descriptor  Riccati  equations  (4.3),  (4.4)  or  their 
time-invariant  counterparts  (4.11),  (4.12).  For  example,  let 


P  = 

^11 

^12 

.  [P"^  +  C'R  ^C]“^  = 

’^11 

1 . 

(N 

rH 

.^12 

^22. 

.''12 
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and  consider  (4.14)  which  in  this  case  reduces  to 


fp,-, 

o' 

'o 

0  1 

11 

.0 

0. 

p 

to 

V  ■  _ _ 

+  BOB’ 


which  is  obviously  inconsistent  with  a  positive  definite  solution  for  P.  Even 
if  one  considers  indefinite  solutions,  we  see  that  none  can  possibly  exist  if 
BOB'  is  not  diagonal.  Indeed  the  only  case  in  which  any  solutions  exist  to 
(4.11),  (4.12)  is  when  BCJB’s  and  C'R  are  both  diagonal.  In  this  case  P  and 
0  are  also  diagonal,  with  the  positive  diagonal  element  corresponding  to  the 
error  covariance  of  the  causal  part  of  the  system  (the  first  state  component) 
and  the  negative  element  to  the  negative  of  the  error  covariance  of  the 
anticausal  (second  state)  component.  Furthermore,  the  diagonal  nature  of  BOB' 
8ind  C’R  implies  that  independent  noises  drive  each  component  and 
independent  observations  are  available  for  each  —  i.e.  the  problem  reduces  to 
the  trivial  and  uninteresting  case  of  two  completely  decoupled  systems. 

On  the  other  hand,  the  generalized  Riccati  equations  (4.14),  (4.15)  admit 
solutions  in  nontrivial  cases.  For  example,  if 


Y 

'l  o' 

B  = 

II 

a 

o 

II 

w 

II 

.1. 

p  1. 

the  solutions  to  (4.14),  (4.15)  are 


'l  o' 

'l  l' 

0  = 

.  ^  = 

.0  3, 

.1  2 
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and 


2  l' 

■3 

o' 

s  = 

,  T  = 

.1  1. 

P 

1, 

This  example  also  illustrates  the  degeneracy  that  arises  in  the  dynamic 
portion  of  the  smoother  for  TPBVDS’s  whenever  either  E  or  A  is  singular. 
Indeed  in  this  case  (4.21)  reduces  to 


5(k+l)  = 


0 

.0 


(5.1) 


-T(k+1)  +  y(k)  (5.2) 

This  is  of  course  2in  extreme  example,  since  the  two  components,  x^  and  Xg,  of 
X  are  essentially  identical  white  noise  sequences  (with  a  sign  inversion  and  a 
one  unit  relative  time  shift)  except  for  the  possible  correlation  between  x(0) 
and  x(K)  introduced  by  the  boundary  conditoins.  However,  while  in  general  the 
system  matrices  in  (4.21)  will  not  be  nilpotent  as  they  are  here,  there  will 
always  be  some  rank  deficiency  if  either  A  or  E  is  singular. 

Filially,  let  us  illustrate  the  rest  of  the  smoothing  solution  for  this 
example.  Even  in  this  degenerate  case  the  one  time-step  delay  between  x^  and 
X2  and  the  nature  of  the  boundary  conditions  can  lead  to  a  nontrivial  form  for 
the  smoother.  In  particular,  suppose  that 


TCk)  = 


0  0 
-1  0 


'l 

o' 

'o 

o' 

1  ! 

'1  1' 

It 

P 

1, 

V  = 

.  _ 

P 

1, 

II 

,1  1. 

(5.3) 
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1 

0 

0 

0 

0 

1 

II 

0 

0 

,  IT  = 

'l  o' 

0 

0 

1 

0 

D 

.0  1. 

p 

0 

p 

1, 

The  dynamics  plus  boundary  conditions  in  this  case  are 

Xj(k+1)  =  u(k) 

k  =  1 . K  -  1 

X2(k)  =  -u(k) 

with  X2(0)  a  unit  variance  random  variable  independent  of  u,  and  with 
X2(K)  =  x^(0)  +  u(0) 

Referring  to  (3.18),  we  have 


f  = 
^11 


0. 

-1 

0 

1 


0 

0 

1 

0 


"V  =0 

,  ^12  -  U 


• 

1 

0 

-1 

-1 

• 

0 

0 

0 

0 

0 

1 

-1 

-1 

IS  - 
'  *22 

0 

1 

0 

0 

0 

0 

0 

0 

1 

0 

1 

0 

p 

0 

0 

1, 

p 

1 

0 

0. 

We  can  then  compute  and  M2  satisfying  (3.19): 


1 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

.  M„  = 

0 

0 

1 

0 

0 

1 

0 

1 

0 

0 

0 

-1 

p 

3 

1 

1. 

-1 

2 

0 

0. 
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The  boundary  conditions  (4.22)  for  (5.1),  (5.2)  then  are 


10  0  0 
0  0  0  0 
0  10  0 


0  0  0  0 

0  0  1  0 

0  0  0  4 

0  0  0  I- 


(5.5) 


where 


“  =  Ml  _i  +  ^2  -1 

K-lIb  i'b  *K’”b  ^b 


Then  applying  (2.16)  to  (5.1),  (5.2),  (5.5)  (with  an  adjustment  for  the 
fact  that  the  smoother  (5.1),  (5.2)  runs  from  1.  rather  than  0,  to  K) ,  we  find 
that 


5(1)  = 


1 

0 

0 

0 

• 

1 

0 

2 

6 

3 

a  + 

2 

0 

11 

0 

11 

11. 

.  11 

5(k)  = 


0  0 
0  1 


y(k-l),  2  <  k  <  K 


0  0 

T(k)  =  y(k)  +  y(k-l).  1  <  k  <  K-2 

-10 


^  0  0  0  o' 

'r(K-l)  =  a  +  y(K-l) 

0-100 
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-rCK)  = 


0 

1 

0 

0 

0 

0 

6 

r  11 

0 

15 

~11 

9 

11. 

a  + 

6 

fll 

0 

nr(l) 


Finally,  using  (4.23)  we  compute 


x(k) 


{Tr(k) 


X(k)  = 


2 

3 

_i 

3 


2  -1 

0  0 


5(k)} 


^(k)  -  t 


2  -1 

-1  2 


k=l.. 

5{k) 


.  .K 


and,  from  (3.22) 


x(0)  = 


5 

8 

1_ 

16 


1_ 

'16 

13 

32 


3_ 

"16 

7_ 

32 


- 

1 

1 

0 

1 

16 

~4 

Ml)  - 

8 

13 

1 

0 

3 

.  32 

8. 

32. 

where 


0 

0 

=  f  ‘ 

11 

...  -1 

+  1f 
^21 

o 

oP 

1 

t _ 

.^K'^b  ^b. 

x(K) 
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Example  5.2:  In  this  example  we  introduce  the  class  of  cyclic  TPBVDS’s  for 
which  the  boundary  condition  (2.2)  takes  the  special  form 
x(0)  =  x(K) 

Equivalently  we  can  think  of  a  cyclic  system  as  being  defined  on  [0,  K-1]  with 
the  boundary  condition 


Ex(0)  -  Ax(K-l)  =  Bu(K-l) 

(so  that  y  in  (2.5a)  is  block-circulant) . 

Consider  the  smoothing  problem  for  such  a  system  when  the  boundary 
measurements  are 

y^  =  Cx(0)  +  r^ 

with  IT,  =  R.  it  is  not  difficult  to  check  that  in  this  case  ^  is  also 
o 

block-circulant  (i.e.  'll =  -s4.  ‘^■^2  ~  "^21  ~  that  the  smoother 

is  also  a  cyclic  TPBVDS  over  [0,  K-1]  (with  no  need  to  move  the  boundary  in 
one  step  as  in  (3.19)  -  (3.22)).  If  we  then  follow  the  procedure  described  in 
Section  4.2,  we  obtain  two  non-descriptor  cyclic  systems 

5(k+l)  =  Fg5(k)  +  Ggy(k)  .  5(0)  =  5(K)  (5.6) 

^(k)  =  F^^(k+1)  +  G^y(k)  .  T(0)  =  t(K)  (5.7) 

where  the  F’s  and  G’s  are  specified  in  (4.21)  and  we  have  adopted  the  notation 
y(0)  =  y(K)  =  y^. 

Obviously  the  symmetry  of  the  cyclic  case  leads  to  some  simplifications. 
In  fact,  note  that  the  two  systems  (5.6),  (5.7),  including  boundary 
conditions,  are  completely  decoupled.  This  greatly  simplifies  their  solution. 


which  we  can  write  as  cyclic  convolutions: 
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K-1 

6(k)  =  I  F/Ggy{k-«-l)  (5.8) 

k=0,l . K-l  • 

K-1 

T(k)  =  [I-F^^]"^  I  F^\y{k+e)  (5.9) 

«=0 

where  we  extend  y(k)  periodically  (i.e.  y(k+K)  =  y(k)).  The  estimate  x(k) 
then  be  computed  from  (4.23a)  again  without  any  need  to  determine  x(0) 
separately  since  we  did  not  need  to  move  the  smoother  boundary. 


can 
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VI.  The  Smoothing  Error  for  TPBVDS’s 

Recall  from  the  development  in  Section  III  that  we  obtained  the  form  of 
the  optimal  smoother  by  expressing  xand  X  in  terms  of  y  and  z  as  in  (3.10) 
and  then  setting  z  to  zero.  Thanks  to  the  orthogonality  of  z  and  y  we  can 
similarly  obtain  an  expression  for  the  smoothing  error  by  setting  z  to  zero  in 
(3.10): 


if  -w' 

•  ■ 

X 

26Qz' 

if' 

mm. 

.  0  . 

(6.1) 


where  x  =  x 


X,  X  =  X  -  X. 


If  we  then  use  these  relationships,  together  with 


(3.7),  (3.8)  we  obtain 


y  -3m'' 

X 

2S 

0 

u 

yv 

.0 

J. 

As  in  Section  III,  this  is  equivalent  to 


(6.2) 


^  ^  r 


'  x(k+l)' 

'  x(k)' 

B  0 

u(k)' 

yN 

=  d 

A 

+ 

—  1 

L-x(k+i)J 

L-x(k)J 

P  C'R  \ 

.r(k). 

(6.3) 


with  boundary  conditions 


x(l) 

-X(l) 


^2*^22^ 


x(N) 
L-X(N) 


Bu(0) 

•  • 

V 

+  M„ 

z 

(6.4) 
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Examining  (6.3),  (6.4).  we  see  that  the  evaluation  of  the  covariance  of 
the  estimation  error  x(k)  corresponds  to  the  computation  of  (the  upper 
left-hand  block  of)  the  covariance  of  the  TPBVDS  (6.3).  (6.4)  driven  by  white 
noise  (u'(k).  r’(k))  and  with  independent  boundary  conditions.  In  the 
Appendix  we  describe  one  method  for  performing  this  computation  for  the 
original  TPBVDS  introduced  in  Section  II.  This  calculation  is  somewhat  more 
complicated  than  the  corresponding  one  for  causal  systems  since  x(k)  in  (2.1) 
is  not  Markov  and  in  fact  is  not  independent  of  future  values  of  u(k) .  We 
refer  the  reader  to  [25]  for  more  on  the  properties  and  calculation  of  the 
covariance  and  correlation  function  of  such  processes. 

We  close  this  section  with  two  final  observations.  First,  note  that  the 
computation  described  in  the  Appendix,  when  applied  to  (6.3),  (6.4)  yields  the 
covariance  of  x(k)  for  k  >  1.  In  order  to  compute  the  covarieince  of  x(0) .  we 
need  to  exeimine  the  counterpart  to  (3.22): 


■  x(0)- 

1 

Bu(0) 

V  1 

.-MO). 

=  D  t 

11 

+  t 

21 

. “\ 

1 

■g 


■x(l)' 

.-Ml). 


"*^21  ”^22^ 


■  x(N)' 

.-MN). 


(6.5) 


The  calculation  of  the  covariance  of  the  left-hand  side  of  (6.5)  then  involves 
the  computation  of  the  covariances  of  and  the  correlations  among  the  various 
random  vectors  appearing  on  the  right-hand  side  of  (6.5).  An  analogous 
computation  is  also  carried  out  in  the  Appendix. 
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The  second  point  concerns  the  diagonal ization  of  (6.2).  In  particular, 
assuming  that  positive  definite  solutions  exist  to  (4.11),  (4.12),  we  can 
perform  analogous  steps  to  those  used  in  Section  4.2  to  transform  (6.2)  into 
the  non-descriptor,  forward  and  backward  pair  of  eqxiations 

^(k)  =  A'S“^Ee"^T(k)  -  A'S"^Bu(k)  +  C‘R"^r(k)  (6.6a) 

6(k+l)  =  AT“^E''l^“^6(k)  +  Bu(k)  +  AT"^‘€*R"^r(k)  (6.6b) 

(with  corresponding,  and  generally  coupled,  boundary  conditions)  with  x(k)  and 

yv  ys.  />.» 

-X(k)  then  obtained  from  (4.23a,  b),  respectively,  with  t  and  5  replaced  by  t 
and  5. 

Equation  (6.6)  is  extremely  useful.  In  the  first  place,  it  provides  the 
forward-backward  decomposition  needed  in  the  covariance  analysis  procedure 
described  in  the  Appendix.  More  importantly,  it  provides  the  basis  for  a 
system- theoretic  investigation  of  the  smoother,  the  initial  parts  of  which  are 
developed  in  the  next  section. 
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VII.  System  Theoretic  Properties  of  the  Smoother 

The  theory  of  filtering  aind  smoothing  for  causal  systems  includes  a  rich 
set  of  system-theoretic  results  related  to  reachability,  observability, 
stability,  eigenstructure,  etc.  Consequently  a  natural  and  important  line  of 
investigation  is  the  development  of  a  parallel  theory  for  TPBVDS’s.  While  no 
such  complete  theory  is  available,  we  can  provide  an  encouraging  start. 
Consider  the  descriptor  system 

Ex(k+1)  =  Ax(k)  +  Bu(k)  (7.1) 

y(k)  =  Oc{k)  (7.2) 

There  are  a  variety  of  notions  eind  definitions  of  reachability  and 
observability  in  the  literature  (see,  for  example  [7,  25]),  but  for  our 
purposes  here,  we  employ  the  following  counterparts  to  one  pair  of  definitions 
used  for  causal  systems. 


Def inition  7.1:  The  system  (7.1)  is  completely  reachable  if  [sE-tAIB]  has 
full  rank  n  for  (s,t)  =  (0,  0).  The  system  (7.1),  (7.2)  is  completely 
observable,  if 

sE  -  tA 
C 

has  full  rank  n  for  (s,t)  =  (0,0). 

Note  that  the  conditions  for  controllability  and  observability  need  only 

9 

be  checked  for  pairs  (s.t)  that  are  ei genmode s  of  the  system,  i.e.  for  which 
det(sE-tA)  =  0. 

9 

We  use  this  definition  of  eigenmodes  as  it  allows  us  to  capture  "eigenmodes 
at  infinity"  (corresponding  to  a  pair  (s,  0))  without  analytic  difficulty. 
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Proposition  7.1-  The  smoothing  error  descriptor  dynamics  (6.3)  are  completely 
reachable  if  and  only  if  (7.1),  (7.2)  is  completely  reachable  and  observable. 
Proof-  Using  the  definitions  of  £  and  in  (3.17)  we  see  that  (6.3)  is 
completley  reachable  if  and  only  if 


sE+tA  -sBQB'  B  0  1 

.-tC'R'^C  -sA’-tE’  0  C'R~\ 

has  full  rank.  Multiplying  on  the  right  by  the  invertible  matrix 


I  0  0  0 

0  I  0  0 

0  sQB’  I  0 

,tC  0  0  1, 


yields 

sE+tA  0  BO' 

0  -sA'-tE'  0  C‘R~^. 


from  which  the  proposition  follows  immediately. 


One  result  that  we  conjecture  is  true  is  that,  as  with  causal  systems  and 
standard  Riccati  equations,  completely  reachibility  and  observability  should 
imply  existence  and  uniqueness  of  positive  definite  solutions  to  the 
generalized  Riccati  equations  (4.14),  (4.15).  One  would  also  expect  that 
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these  conditions  would  imply  filter  stability.  We  can  prove  one  result  along 
these  lines. 

Proposition  7.2:  Suppose  that  positive  definite  solutions  exist  to 
(4.14)  and  (4.15),  and  suppose  also  that  (7.1),  (7.2)  is  completely  reachable 
and  observable.  Then  the  smoother  is  forward-backward  stable. 

Proof '  What  we  wish  to  examine  is  the  stability  of  (6.3)  or, 
equivalently,  (6.6a)  and  (6.6b).  For  the  latter  equations  we  can  write  down 
standard  Lyapunov  stability  equations: 

Pg  -  (AT"^E’'i^“^)Pg(AT“^E''J'"^) '  =  BOB'  +  AT“^C’R“^Cr"^A*  (7.3) 

P^  -  (A'S“^E9"^)P^(A‘S"^E9"^)’  =  AS“^BQB*S~^A  +  C'R"^C  (7.4) 

By  Proposition  7.1,  (6.3)  is  completely  reachable,  and  therefore  (6.6a)  and 
(6.6b)  are  each  completely  reachable.  Therefore  the  forward-backward 
stability  of  (6.6a)  and  (6.6b)  is  equivalent  to  the  existence  of 
positive-definite  solutions  to  (7.3),  (7.4).  However,  examination  of  (4.16) 
and  (4.17)  shows  that  the  solutions  to  (7.3),  (7.4)  are 

Pg  =  ^  P^  =  9  (7.5) 

which  yields  the  result. 

This  result  deserves  some  comment.  First,  recall  from  Section  II  that 
the  construction  of  the  Kronecker  canonical  form  is  one  general  method  for 
constructing  a  forward-backward  stable  decomposition  of  a  general  TPBVDS. 

What  Proposition  7.2  describes  is  a  second-way  in  which  to  accomplish  this  for 
Hamiltonian  TPBVDS.  Second,  given  what  we  know  about  the  causal  case,  it  is 
not  surprising  that  there  is  a  close  connection  between  generalized  Riccati 
equations  and  Hamiltonian  eigenstructure.  Indeed  one  might  expect  there  to  be 
a  generalized  Hamiltonian  eigenvector  approach  to  solving  these  equations  that 
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is  analogous  to  the  popular  method  for  standard  Riccati  equations  [27.  28], 
Such  a  development  remains  for  the  future,  but  we  can  derive  a  related  result: 

Proposition  7.3:  If  (Sq,  t^)  is  an  eigenmode  of  the  pencil  {S ,  d,}  then 
so  is  (tQ,  Sq) 

Proof :  Note  first  that  if  (Sq,  0)  is  an  eigenmode,  i.e.  if 


det(«)  = 


-BCB’’ 
-A'  . 


=  0 


then  (0,  Sq)  is  also  an  eigenmode,  i.e. 


det(s4) 


det 


A 

C‘R"^C 


0 

E’ 


=  0 


Consider  then  any  eigenmode  (Sq,  tQ)  with  Sq,  ^q^O.  The  following 
computation  then  shows  that  (tQ,  Sq)  is  also  an  eigenmode: 
det(tQS-SQj4)  =  det(tQS ’-SQsi' ) 


0 

1/t  I 

0 

=  det 

U 

(tQr-SQ^') 

u 

I/SqI 

0 

rt 

O 

0 

=  det  (sQ^-tQsi)  =  0 

Note  that  this  is  the  generalization  of  the  usual  reciprocal  symmetry  of 
Hamiltonian  eigenvalues.  From  this  result  we  can  immediately  deduce  that  the 
system  matrice  AT  ^  and  A'S  ^E0  ^  associated  with  the  forward-backward 

smoother  decomposition  have  identical  eigenvalues. 
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Finally  we  present  one  additional  result. 


Proposition  7.4:  The  pencil  (5,  sd)  is  forward-backward  stable  if  {E,  A)  is. 
Proof :  We  need  to  show  that  (l.e"^^)  is  not  an  eigenmode  of  (S,  sd)  for  any  u. 
Consider 


'  E-e-^’^A  -BOB' 

-A’+e^^^E', 


(7.6) 


Since  (E,  A)  is  forward-backward  stable 
r  =  E  -  eJ“A 

is  invertible.  Therefore  the  invertiblity  of  (7.6)  follows  if  we  can  show 
that 

+  C*R"^Cr"^BQB' 

(where  "H"  denotes  conjugate  transpose),  or  equivalently  I  +  MK  is  invertible, 
where 

M  =  C'R"^C.  K  =  (r"^B)Q(r"^B)^ 

This  follows  from  the  positive  semi-definiteness  of  M  and  K  and  the  consequent 
nonnegativity  of  the  eigenvalues  of  MK. 


This  result  roughly  corresponds  to  the  causal  result  stating  that  the 
Kalman  filter  is  stable  if  the  original  system  is,  independent  of  any 
controllability  eind  observability  results.  What  we  conjecture  is  also  true  is 
a  blending  of  Propositions  7.2  and  7.4,  namely  that  the  smoother  is 
forward-backward  stable  if  the  system  (7.1),  (7.2)  is  forward-backward 
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stabilizable  and  detectable,  i.e.  if  [sE-tA!B]  and  [sE'-tA'iC]  have  full  rank 
for  all  eigenmodes  such  that  |s/t|  =  1. 
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VIII.  Conclusions 

In  this  paper  we  have  investigated  the  optimal  estimation  problem  for 
two-point  boundary-value  descriptor  systems  (TPBVDS’s).  Using  the  method  of 
complementary  processes  we  developed  a  generalization  of  the  Hamiltonian  form 
of  the  optimal  smoother  for  causal  systems.  This  genrealized  Hamiltonian 
system  is  itself  a  TPBVDS.  In  addition,  we  have  generalized  the  notion  of 
Hamiltonian  diagonal ization  as  a  method  for  reducing  the  smoother  to  two 
systems  of  lower  order.  Both  of  the  approaches  described  involve 
generalizations  of  standard  Riccati  equations.  One  of  these,  corresponding  to 
descriptor  Riccati  equations  that  have  appeared  in  the  literature,  is  shown  to 
work  only  in  certain  cases  and  is  not  appropriate  when  the  system  dynamics  are 
intrinsically  acausal,  i.e.  when  both  system  matrices  E  and  A  are  singular. 
However,  our  second  approach,  involving  what  we  call  generalized  Riccati 
equations,  appears  to  offer  much  promise.  Indeed  we  have  illustrated  that  it 
does  provide  a  viable  approach  in  the  acausal  case.  Furthermore,  the  results 
presented  in  Section  VII  indicate  that  there  is  likely  to  be  a  complete  system 
theory  for  these  new  Riccati-like  equations  end  the  associated  generalized 
Hamiltonian  system. 

There  are  numerous  open  questions  raised  by  the  work  described  in  this 
paper.  In  the  previous  section  we  indicated  several  of  these,  namely 
existence  and  uniqueness  conditions  for  the  generalized  Riccati  equations, 
Hamiltonian  eigenvector  solutions  to  these  equations,  and  weaker  stability 
conditions  involving  stabilizability  and  detectability.  Also,  an  important 
question  is  the  relationship  of  the  solutions  of  these  Riccati  equations  to 
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the  estimation  error  covariance,  whose  computation  we  can  presently  describe 
only  in  the  mechanical  manner  given  in  Appendix  B.  In  addition,  there  are 
other  important  questions  related  to  alternate  notions  of  stability  and  weaker 
concepts  of  controllability  and  observability  that  make  sense  for  TPBVDS’s  and 
that  are  developed  in  [25].  For  example,  consider  the  cyclic  system  described 
in  Section  VI.  Such  a  system  can  be  thought  of  naturally  as  living  on  a 
discretized  version  of  the  circle.  Forward-backward  stability  in  this  case 
corresponds  to  clockwise  and  counter-clockwise  staiblity.  An  alternate  notion 
of  recursion  developed  in  [25]  involves  computations  that  begin  at  one  point 
and  proceed  simultaneously  in  clockwise  eind  counterclockwise  directions  until 
the  entire  circle  is  covered.  In  this  case  stability  would  correspond  to 
convergent  behavior  as  the  radius  of  the  circle  grows  without  bound.  As 
described  in  [25]  it  is  possible  to  develop  a  stability  theory  and  in  fact 
generalized  Lyapunov  methods  along  these  lines.  The  implications  of  these 
concepts  for  the  smoother  represents  another  intriguing  line  of  investigation. 
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Appendix:  Covariance  Analysis  for  TPBVDS’s 

In  this  appendix  we  develop  formulas  for  covariance  analysis  of  TPBVIS’s. 
As  a  starting  point  for  this  computation,  we  assume  that  our  TPBVDS  has  been 
placed  in  the  forward-backward  form  given  in  (2.9),  (2.11).  The  general 
solution  for  this  system  is  given  in  (2.16).  Given  the  independence  of  the 
boundary  value  v  and  the  white  sequence  u(k) ,  we  see  that  the  covariance  of 
x(k)  can  be  expressed  in  terms  of  the  covariance,  U^,  of  v  and  the  three 


quantities 

P^°(k)  =  ECXj°(k)x^°(k)’]  (B.la) 

P^°(k)  =  E[Xj^°(k)x^°(k)-3  (B.lb) 

P^^°(n.k)  =  E[x^°(n)Xj^°(k)-]  (B.lc) 

The  computations  of  these  quantities  are  straightforward: 

P^°(k+1)  =  A^P^°(k)A^’  +  B^QB^'  Pf°(0)  =0  (B.2a) 

P^°(k-1)  =  Aj^P^°(k)A^^*  +  B^QB^'  P^°(K)  =  0  (B.2b) 

and 


Pj.^°{n,k)  =  \ 


,  n  <  k 


IT°fb(n)(Ab’)''  ^  -  Af''“\b°(k)  .  n>k 


(B.2c) 


where 

ir^^°(k+l)  =  A^irj^°(k)Aj^'  +  B^QB^' 

Given  these  quantities,  we  can  now  determine  an  expression  for 


x^(k) 


2(k)  =  E  { 


[x^(k)J 


[x^*(k),  x^‘(k)]} 


(B.3) 
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2(k)  =  G(k)IT^G‘(k)  +  G(k)\J^(k)  +  '^'•(k)G•(k)  +  G(k)AG'(k) 


0 


0  P^°{k) 


(B.4) 


where 


G(k) 


(B.5) 


^(k)  =  -  ^p^°(k);p^^°(K.k)]  +  v^_Q[P^^°(k.0)-:A^\°(k)] 

(B.6) 


A 


[V 


■  V 

f.K.  b.O- 


Pf°(K) 

Pft,°(K.O) 


Pfb  (K.o) 


fv  *1 
^f  ,K 

V  ■ 

J 

L  b.oj 

{B.7) 


As  mentioned  in  Section  6.  the  computation  of  the  error  covariance  at  the 
initial  point  in  the  interval  of  interest  involves  an  additional  computation. 
In  the  remainder  of  this  Appendix  we  describe  the  corresponding  calculation 
for  (2.9).  (2.11).  Specifically,  suppose  we  would  like  to  compute  the 
covariance  of 


Xf(0) 

Xj^(O) 


N, 


Xf(K) 

x^(K) 


T3  = 


(B.8) 
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where  f  is  a  zero-mean  random  vector  correlated  with  the  boundary  condition  v 
but  independent  of  u.  Let 

E[ff']  =  P^  ECfv’]  =  E[T7n']=P^  (B.9) 

Then,  with  the  help  of  (2.16)  we  have 


P  =  N,PeN, ' 

+ 

Pf(0) 

0 

No' 

+  No 

“  ■ 

Pf(K)  0 

T?  1  f  1 

2 

0 

» 

P^CO) 

2 

3 

0  P^(K) 

+  N^P^^G-(0)N2'  +  N^GCOP^^’Nj-  H-  NjP^^G'(K)N3' 


+  N3G(K)P^^’N^*  +  N2P’(K.0)N3'  +  N3P(K.0)N2' 


where 


P(n,k)  =  E{ 


’Xf(n)' 


Cx^’(k)  x^'{k)]} 


,  n  >  k 


(B.IO) 


(B.ll) 


can  be  calculated  in  the  same  manner  as  2(k): 


P(n.k)  =  G(n)I7^G'(k)  +  G(n)'li-(k)  +  'l>’(n)G'(k) 


+  G(n)AG’(k)  + 


0 


(k) 


0 


0 


(B.12) 
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